import random
import cv2
import numpy as np
import matplotlib.pyplot as plt
plt.rcParams['axes.unicode_minus']=False
plt.rcParams['font.sans-serif']=['SimHei']
import time
from PIL import Image
import torchvision.transforms as transforms
import os


m=2    #隶属度函数因子
his_clust_num=[]

class FCM:
    def __init__(self, data, clust_num,iter_num):
        self.data = data
        self.cnum = clust_num  #聚类中心数
        self.sample_num=data.shape[0]  #样本集
        self.dim = data.shape[-1]  # 数据最后一维度数  特征数
        Jlist=[]   # 存储损失函数计算值的矩阵

        U = self.Initial_U(self.sample_num, self.cnum)  #初始化隶属度矩阵
        for i in range(0, iter_num):
            C = self.Cen_Iter(self.data, U, self.cnum)
            U = self.U_Iter(U, C)
            if (i==iter_num-1):
                print("第%d次迭代" %(i+1))
                print("聚类中心",C)
            J = self.J_calcu(self.data, U, C)  # 计算损失函数
            Jlist = np.append(Jlist, J)
        self.label = np.argmax(U, axis=0)  # 所有样本的分类标签
        self.Clast = C    # 最后的类中心矩阵
        self.Jlist = Jlist  # 存储损失函数的矩阵

    # 初始化隶属度矩阵U
    def Initial_U(self, sample_num, cluster_n):
        U = np.random.rand(sample_num, cluster_n)  # sample_num为样本个数, cluster_n为分类数
        row_sum = np.sum(U, axis=1)  # 按行求和 row_sum: sample_num*1
        row_sum = 1 / row_sum    # 该矩阵每个数取倒数
        U = np.multiply(U.T, row_sum)  # 确保U的每列和为1 (cluster_n*sample_num).*(sample_num*1)
        return U   # cluster_n*sample_num

    # 计算类中心
    def Cen_Iter(self, data, U, cluster_n):
        c_new = np.empty(shape=[0, self.dim])  # self.dim为样本矩阵的最后一维度
        for i in range(0, cluster_n):          # 图片像素值的dim为1
            u_ij_m = U[i, :] ** m  # (sample_num,)
            sum_u = np.sum(u_ij_m)
            ux = np.dot(u_ij_m, data)  # (dim,)
            ux = np.reshape(ux, (1, self.dim))  # (1,dim)
            c_new = np.append(c_new, ux / sum_u, axis=0)   # 按列的方向添加类中心到类中心矩阵
        return c_new  # cluster_num*dim

    # 隶属度矩阵迭代
    def U_Iter(self, U, c):
        for i in range(0, self.cnum):
            for j in range(0, self.sample_num):
                sum = 0
                for k in range(0, self.cnum):
                    temp = (np.linalg.norm(self.data[j, :] - c[i, :]) / np.linalg.norm(self.data[j, :] - c[k, :])) ** (2 / (m - 1))
                    sum = temp + sum
                U[i, j] = 1 / sum
        return U

    # 计算损失函数值
    def J_calcu(self, data, U, c):
        temp1 = np.zeros(U.shape)
        for i in range(0, U.shape[0]):
            for j in range(0, U.shape[1]):
                temp1[i, j] = (np.linalg.norm(data[j, :] - c[i, :])) ** 2 * U[i, j] ** m

        J = np.sum(np.sum(temp1))
        print("损失函数值:%.2f" %J)
        return J

def FCM_pic_cut(img_path,gray=False,clustercenternum=3,iternum=5, outpath=None):
    if gray:
        img=cv2.imread(img_path,0)  #灰度图
        data=img.reshape(img.shape[0]*img.shape[1],1) #将图片拉成一列

    else:
        img=cv2.imread(img_path)
        img=cv2.cvtColor(img,cv2.COLOR_BGR2RGB) #转化为RGB，不然plt时图片会怪怪的
        data=img.reshape(-1,3)  # 将三维降成二维

    print("开始聚类")
    start = time.time()
    test=FCM(data,clustercenternum,iternum)
    cluster=test.label  # 聚类结果
    center=test.Clast # 聚类中心
    his_clust_num.append(test.Jlist[-1])
    end = time.time()

    print("聚类完成，开始生成图片")
    #print(cluster)
    new_img=center[cluster] # 根据聚类结果和聚类中心构建新图像
    new_img=np.reshape(new_img,img.shape) #矩阵转成原来图片的形状
    #print(new_img)
    new_img=new_img.astype('uint8')  # 要变成图像得数据得转换成uint8
    plt.show()
    new_img = Image.fromarray(np.uint8(new_img))
    out_img = new_img.convert('L')
    print('FPS:',1/(end - start))
    plt.imsave(outpath,out_img) # 保存图片

if __name__ == '__main__':
    image_path = r"E:\Image2\G50F008021_2_0_1_0.tif"
    save_dir = r"E:\Yang\image"
    file_name = f"result.jpg"
    save_path = os.path.join(save_dir, file_name)
    FCM_pic_cut(image_path,gray=False,clustercenternum=3,outpath=save_path)

    # image_path = r"E:\Image2\G50F008021_2_0_1_0.tif"
    # save_dir = r"E:\Yang\image"
    # # image_path = "190811vvaoligti.jpg"
    # # 读取图片
    # image = Image.open(image_path)
    # # 输出维度
    # print("RGB图像的维度：", np.array(image).shape)
    # # 显示原图
    # image.show()
    # # RGB转换灰度图像
    # image_transforms = transforms.Compose([transforms.Grayscale(1)])
    # image = image_transforms(image)
    # # 输出灰度图像的维度
    # print("灰度图像维度： ", np.array(image).shape)
    # # 显示灰度图像
    # file_name = f"result.jpg"
    # save_path = os.path.join(save_dir, file_name)
    # print(save_path)
    # image.save(save_path)

